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Abstract 

The effects of confinement on colloidal self-assembly in the case of hxed number of confined 
particles are studied in the one dimensional lattice model solved exactly in the Grand Canonical 
Ensemble (GCE) in [J. P§kalski et ah J. Chem. Phys. 142, 014903 (2015)]. The model con¬ 
siders a pair interaction defined by a short-range attraction plus a longer-range repulsion. We 
consider thermodynamic states corresponding to self-assembly into clusters. Both, fixed and adap¬ 
tive boundaries are studied. Eor fixed boundaries, there are particular states in which, for equal 
average densities, the number of clusters in the GCE is larger than in the Canonical Ensemble. 
The dependence of pressure on density has a different form when the system size changes with fixed 
number of particles and when the number of particles changes with fixed size of the system. In the 
former case the pressure has a nonmonotonic dependence on the system size. The anomalous in¬ 
crease of pressure for expanding system is accompanied by formation of a larger number of smaller 
clusters. In the case of elastic confining surfaces we observe a bistability, i.e. two significantly 
different system sizes occur with almost the same probability. The mechanism of the bistability 
in the closed system is different to that of the case of permeable walls, where the two equilibrium 
system sizes correspond to a different number of particles. 
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I. INTRODUCTION 


Confinement has a significant effect on fluids when the separation between the confining 
surfaces is comparable with the characteristic structural length of the confined fluid jJ, l2|. 
In simple fluids such pore thicknesses are of order of Angstroms. The effects of confinement 
depend on whether the boundaries are rigid or adaptive. In the latter case mechanical 
equilibrium between the stress in the solid walls resulting from the swelling or shrinking of 


the pore, and the solvation force 


,3 


induced on the walls by the confined fluid can lead 


to contraction or expansion of the pore 3- Similar dependence on the type of confinement 
is expected in the case of self-assembling systems, but on different length and energy scales 


in| . The packing effects of molecules are replaced in this case by the packing effects of 
micelles, clusters or layers that are much larger and softer. 

In this work we focus on colloidal self-assembly in thermodynamic conditions correspond¬ 
ing to self-assembly into clusters or layers in the bulkjs, 11-15|. We study rigid and elastic 
boundaries, and consider separations of the confining walls comparable with a few struc¬ 
tural units. The adaptive confinement in this case means confining surfaces that are soft 
and separated by hundreds of nanometers. Important examples of such a confinement exist 
in organella or in lipid or polymeric vesicles. One can expect that the shape of the outer 
membrane and the ordering inside the vesicle or organellum influence each other in a way 
that depends on the elasticity of the membrane. We shall consider self-assembly in systems 
with fixed boundaries and with boundaries with different elasticity, from stiff to very soft. 

Some of the membranes are permeable, while some other are not. Thus, a question arises 
if fluctuations of the number of confined particles have any effect on the properties of a 
self-assembling system confined by rigid or adaptive boundaries. This question motivates 
our present study. We ask how the effects of confinement on the self-assembling system 
depend on the contact with a reservoir of particles. 

The above problems are very difficult for realistic models of self-assembling systems con¬ 
fined by elastic boundaries. In order to gain some first insight, however, one can consider 
simplified models. In this work we consider the simple generic model of the system with 


competing short-range attraction and longer-range repulsion (SALR) that can be so 


actly 


16| . The model was solved exactly and systematically analyzed in the bulk 


16 


ved ex¬ 


in confinement 


181 and 


19| in the grand canonical ensemble (GCE). Recently a procedure to extract 
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full Canonical information from Grand-Canonical results has been proposed in Ref. 20 |. 
In principle, it is possible to apply this technique to our exact GCE results. However, in 
the current problem the presence of energetic terms and low temperature pose numerical 
difficulties in the mapping procedure between the results in both ensembles. For this reason 
we decide to perform MC simulations in the canonical ensemble (CE). 

In the case of hxed boundaries we address the question of how the fluctuations of the 
positions and sizes of the clusters (their dynamical assembly and dissociation) are coupled 
with the fluctuations of the number of particles in the system. In the case of inhomogeneous 
distribution of particles it is not a priori obvious that the largest fluctuations in the total 
number of particles lead to the largest differences between the density prohles in the canonical 
and grand canonical ensembles. We shall compare the density prohles and the pressure in 
the GCE and CE ensembles, with the average number of particles in the GCE equal to the 
number of particles Nq in the CE. We shall pay particular attention to values of Nq that 
are too small or too large for a given system size L for formation of periodically distributed 
layers of particles that minimizes the system energy. Roughly speaking, in the SALR systems 
the minimum energy is assumed when the individual clusters have the lowest energy (no 
intra-cluster repulsion), and are separated by the smallest distances corresponding to no 
inter-cluster repulsion. When the number of particles is too small or too large for a given L 
for formation of the optimal bulk structure, some structural deformations must occur. Our 
purpose is to compare these deformations and their effect on the mechanical properties in 
the GCE and CE. 

Adaptive boundaries were studied in Ref. 2l| in the case of hxed number of discs sur¬ 


rounded by particles kept by laser tweezers. The system exhibited bistability; either hexag¬ 
onal structure (and modihed boundaries) or concentric rings of particles occurred. We 
observed a different bistability in Ref. 19|, where we studied particles interacting with the 
SALR potential and conhned by boundaries, whose separation could be varied at the cost of 
elastic energy. When the average number of particles is too small or too large for given L for 
formation of the optimal bulk structure, some compromise must occur between the excess of 
the free energy associated with the structural deformations, and the elastic energy cost due 
to adjusting the size of the compartment to the optimal structure of the conhned huid. We 
obtained exact results for a one dimensional (Id) model in GCE, and required mechanical 
equilibrium between the huid-induced solvation force and the elastic force present when the 
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system shrinks or expands. We fonnd that when the equilibrium width of the empty con¬ 
tainer corresponded to the largest stress of the conhned colloidal system, two signihcantly 
different system sizes were almost equally probable. The large size fluctuations are accom¬ 
panied by absorption or evaporation of a whole cluster. Clearly, for hxed number of particles 
such large size fluctuations are not possible. Hence, the permeability plays an important 
role in the case of elastic boundaries. In this work we verify if in the case of hxed number 
of particles the bistability of the system size can still exist due to some other mechanism. 

We present the model and the simulation methods in Sec. HIl The density prohles are 
described and compared with the results obtained in the GCE in Sec. IIIIl The mechanical 
properties are discussed in Sec. HVl We compare the dependence of pressure on density in 
the CE (hxed N and varying L) and in the GCE for several hxed values of L and varying 
(TV). The dependence of pressure on the system size for hxed N or (N) in the CE or GCE 
respectively is also discussed. In Sec. |V] we consider elastic boundaries and compute the 
average system size as a function of N for various elastic constants of the walls. For selected 
cases the histograms for the huctuating width of the system are presented. Sec. IVII contains 
a short summary and concluding remarks. 


II. THE MODEL AND THE METHODS 


In this section we briehy describe the model and the methods used for its analysis. More 
detailed descriptions of the model and the transfer matrix method used for hnding the exact 


solutions of the model in the GCE can be found in 


19| 


A. The model 

We consider a one-dimensional (Id) lattice model and assume an isotropic ehective pair 
interaction potential of the SALR type. The ranges of the competing attractive and repulsive 
parts are chosen such that small clusters are formed for some range of thermodynamic 
variables. Namely, we assume that the nearest neighbors attract each other and the third 
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neighbors repel each other, i.e. the effective pair interaction potential is 


V{Ax) = < 


— Ji for |Aa;| = 1, 
+J 2 for |Aa:| = 3, 
0 otherwise, 


( 1 ) 


where Ji is the energy of attraction, J 2 is the energy of repulsion and the unit distance is the 
particle diameter a, and Ax is the distance between two sites of the lattice. Such a shape of 
the effective potential can be found for charged particles in solvents that induce short range 


attraction between the particles, e.g. for lysozyme molecules in water 
a solvent containing short-chain polymers 12j. 


22l | or for colloids in 


In Ref. 


16l-ll8| we assumed periodic boundary conditions and extensively investigated the 


spontaneous pattern formation of particles in the bulk. In order to study the effects of 
confinement on the self-assembled structures, we changed the boundary conditions of the 
model from periodic to rigid or elastic in Ref. [^. In the current study we continue the 
investigation of the confined system with particles interacting via the pair potential given 
by Eq. ([T]). The confinement is assumed to be electrically neutral, hence the interaction of 
the particles with the walls is short range. The Hamiltonians in the canonical, and Grand 
canonical ensembles take the forms: 

, L L 


p{x)V{x — x')p{x') + hip{l) + hLp{L) {iixed iV); 


x=l x'=l 


H[{p}] = U[{p}] - p^p{x). {iixed p), 


( 2 ) 


( 3 ) 


X=1 


where hi and are the energies of the interactions between the confining walls and the 
particles located at the 1-st and the L-th site respectively, p{x) is the occupancy operator 
such that p{x) = 1 if the x-th site is occupied and p{x) = 0 otherwise. The number of 
particles N (fixed in the canonical ensemble) is given by: 




( 4 ) 


X=1 


As in Ref. 


19| . we assume that the confinement can be either rigid or elastic. For rigid 


boundary conditions the distance L between the confining walls is fixed, while in the case 
of elastic walls we assume that it may oscillate around Lq, which is the equilibrium width 
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of an empty system. The Hamiltonian of the system with elastic boundary conditions is for 
the case of permeable walls: 

n = H + k{L- LoY, (5) 

where k is the elastic constant, and for the fixed number of particles the Hamiltonian takes 
the form: 

U = U + k{L- Lof. (6) 

We set the energy of attraction Ji as energy unit and introduce the following dimensionless 
variables 

T* = ksT/J,, (3* = hlknT, J* = J 2 /J 1 , (7) 

hi = hi/Ji, hl = hL/Ji, /i* = /i/Ji (8) 

where ks is the Boltzmann’s constant and T is the temperature. From now on we set J* = 3, 
and hi = h*^ = ±1, and study the role of the temperature and the chemical potential or the 
number of particles. For the selected parameters the energy ([2]) takes the minimum when 
clusters composed of 3 particles are separated by 3 empty sites, and a cluster is attached to 
each boundary. Such an optimal structure is possible only for L = 2N — 3. 


B. The computational method 


In order to hnd the exact solution of the model in the GCE for T > 0 the transfer matrix 
method was applied. We have found exact expressions for the partition function and the 
local density. The details of the derivations are provided in the the appendix of Ref. 19 |. 
Here we will present only the final formula for the partition function with a brief description 
of the notation. The partition function for the system of size L = 3M + j, with M a natural 
number and j = 0,1, 2 is given by 


5(1) S(M) 


(9) 


where S{n) = {p{3n — 2), p(3n — 1), p{3n)) with n = 1,... M and T is the transfer matrix 
with the matrix elements defined as 


T(^(n), S{n + 1)) = e-^*^R^G).s(n+i)]^ 


( 10 ) 
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where 


3n 

H*[S{n),S{n + l)]= ^ [ - p{x)p{x + 1) + J*p{x)p{x + 3) - fi*p{x)]. (11) 

x=3n—2 


depend on j and the rather lengthy expressions are provided in Appendix 1X1 
Having a formula for the partition function one can derive exact expressions for the pressure 
and the local density. For more details as well as for asymptotic expressions of the exact 


solutions see 


19|. 


C. The simulation methods 

In the current study we compare the exact solution obtained in GCE with results of 
MC simulation in CE. In addition to the simulations in the CE we also carried out some 
runs in the GCE with the aim of cross-checking the consistency between the simulation 
codes and the numerical treatments based on the transfer matrix methods. The simulation 
procedures make use of the Metropolis criterion implemented for two kinds of MC steps: 
translations of one particle in the Canonical Ensemble and particle insertions or deletions 
in the Grand Canonical ensemble. For each step a trial conhguration is generated and it is 
accepted with probability: min [1, exp{—f3AH)], where AH is the change of the Hamiltonian 
in the trial step. The trial conhguration for a particle translation is generated by moving 
a randomly chosen particle to a randomly chosen empty site, which is equivalent to choose 
two sites of the lattice x and x' so that p{x) ^ p{x'), and interchange their occupancy states, 
so that p^^^°‘\x) = p{x'), and /}*''*“* (x') = p{x). A Monte Carlo step in the Grand Canonical 
simulations implies the insertion or deletion of one particle, this is achieved by choosing 
at random one site of the lattice, x, and generating the trial conhguration by flipping its 
occupancy state from its current value p{x) to the trial value (x) = 1 — p{x). 

The computation of pressure from simulation of lattice models is usually carried out by 
means of the integration of the Grand Potential in the GCE, because the known relation for 
the canonical ensemble: flp = — {d{flA)/dV)^j, (where A is the Helmholtz free energy, and 
V the generalized volume) is hard to translate into an efficient numerical procedure due to 
the discreteness of the volume in lattice systems. We have found, however, that for our Id 
system it is feasible to compute the pressure p{N, L, T) in the CE, by an algorithm based 
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on the discretization of the derivative of A{N, L,T) with respect to the system size as: 


p±{N, L, T) = T[A(iV, L ± 1, T) - A{N, L, T)]. (12) 


The pressure, in the terms of the canonical partition function Q{N, L,T) can be written as 


P±{N,L,T) 


1 Q{N,L±1,T) 

P "" Q{N,L,T) 


(13) 


The two ways of discretization, p_|_ and p_, lead to two different methods of computing 
the pressure, the virtual expansion and the virtual contraction respectively, from a direct 
analysis of the conhgurations from a simulation run at conditions {N,L,T). 

In the virtual expansion an empty site is added at a randomly chosen position of the 
system. For a conhned system of size L, there are L + 1 possibilities of performing such 
an insertion, namely L — 1 cases where the inserted site is located between two sites of the 
system plus two insertions between the walls and the hrst or last site. Considering the L + 1 
possible ways of inserting an empty site on each of the microstates of the system of size L 
with N occupied sites (with 0 < iV < L), we get L + 1 — iV identical copies of each of the 
the microstates of the system of size L + 1 and N occupied sites (see Appendix [B]) . Let us 
denote by pi a given conhguration of the system with N particles and L sites, with potential 
energy given by U{pl). If we dehne pL+i{pL,k) as the conhguration with N particles of a 
system with L + 1 sites built from ph by inserting a site at position fc, and denote by 
the sum over all possible microstates of this system with N particles and L sites, then we 
can write Eq. (Hi as: 


p+{N,L,T) 


1 exp [-PU{pl)] ELo^xp [-pU{pL+i{pL, k)) + I3 U{pl)] 

P ^ (L + 1-A^)^-^exp[-/317(pL)] 


= ksT 111 ^ ^ ^ ^ exp [-^AUins] ^ , (15) 

where {■)l is the average value of ■ when sampled on a system of size L, and AUins is the 
difference between the energies of the systems with L + 1 and L sites. Analogously, a formula 
for pressure p- computed by the virtual contraction scheme can be derived 


p_ {N, L,T) = -ksT In ^ ^ ^ ^ exp [-/3AUdei] ^ , (16) 

where AUdei = f^(PL-i) — U{pl) is the variation of energy when a conhguration /+_i, of 
N particles and L — 1 sites is generated by eliminating one of the empty sites from the 






configuration px of a system with N particles and L sites. The virtual contraction method 
is inefficient at high densities, therefore we used it only for verihcation of the results obtained 
via the virtual expansion method, since by construction we expect: 

p^{N,L.T)=p_iN,L + l,T). (17) 


In order to calculate properties of the system with elastic boundary conditions described 
by the Hamiltonian U given in Eq. ([H]), one needs to perform two additional types of MC 
steps. The hrst one is the intercalation of an empty site into a randomly chosen place of the 
system, the second one is the removal of a randomly chosen empty site of the system. The 
acceptance probability of the hrst move for system with L sites and N particles is 

A{L + 1\L) = min|l,exp[-/jA7/] ^ (18) 


where AU is the change of the energy after the size modihcation. Analogously the probability 
of acceptance of the move in which an empty site is removed is given by: 


A{L-1\L) 



(19) 


Notice that this procedure resembles a lattice version of isothermal-isobaric {NpT) simu¬ 
lation, in which the energy term introduced through the elastic force plays the role of the 
external held. 


III. DISTRIBUTION OF PARTICLES BETWEEN RIGID WALLS 

In this section we consider a system containing Nq particles between rigid walls separated 
by a hxed distance L. The main question is how the particles self-assemble if Nq is such 
that the equilibrium bulk structure is not possible. The distribution of particles for hxed Nq 
will be compared with the distribution of particles in the open system, where the number 
of particles N huctuates in such a way that (N) = Nq. 

In Fig{T] the density prohles obtained by the MC simulations in the CE are compared 
with the exact results obtained in the GCE by the transfer matrix method described in 
Rel,Q, We cw . ^ 50, T* ^ 0,3 eeve.a, values ot 1V„, In each case IKe chennca, 
potential in the GCE was hxed to the value that corresponds to {N) = Nq. We used the 
exact expression for density as a function of the chemical potential that was obtained in 
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Ref.[l^. For L = 51 the optimal number of clusters for the considered range of fj* is 9 
(hence N = 27), since the sequence of three occupied sites followed by three empty sites 
can be formed, with two clusters adsorbed at the attractive surfaces. The energy for such 
a structure assumes a minimum (there are as many attracting pairs as possible with no 
repulsion). For L = 50 only a small defect in the ordered structure occurs. We can see a 
very good agreement between the two ensembles for small as well as for large number of 
particles. For iV = 24 corresponding to eight clusters, however, the number of maxima in 
the GCE is larger than in the CE. This result is even more surprising when we consider 
the fluctuation of the number of particles in the GCE (Fig|2]). One can see that the largest 
discrepancy between the density profiles does not occur for the largest fluctuation of the 
number of confined particles in the GCE. 

In order to understand why the distributions of the particles in the two ensembles are 
different when one cluster in the CE “is missing”, let us consider the ground state (GS), 
T* = 0. The microstates present in the GS for L = 15 and iV = 5, 6, 7 are shown in FigJSl 
The GS in Figj3] shows that even small fluctuations of the number of particles - addition of 
one particle in our case - can lead to a significant change in the distribution of the particles. 
This is the case when (iV) is a multiple of 3, and there is a free space for an extra cluster 
(with no cluster-cluster repulsion). When one additional particle enters the system, the 
energy change is J* — 1 when one of the clusters grows to contain 4 particles, or 0, when the 
new particle is sufficiently far from the clusters, or one of the clusters together with the new 
particle form two clusters composed of two particles. In FigJSl bottom row, the two latter 
cases are shown. The states shown in Fig|3]are energetically favorable for J* > 1. Thus, in 
the GGE with (iV) = 6 such microstates will appear quite often. As a result, an additional 
maximum in the average density profile occurs. 

The above simple considerations show that spatial distribution of particles in the GE and 
GGE can be qualitatively different. This qualitative difference is not present for the largest 
fluctuation of the number of particles in the GGE. Even a small fluctuation of the number 
of particles can lead to a change of the number of clusters, because the sizes of the clusters 
can fluctuate, especially when N is not a multiple of 3. When {N) is small enough, the 
additional clusters can occupy the empty space and no inter-cluster repulsion appears. 
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FIG. 1: Comparison of the GCE (dashed line) and the CE (solid line) density profiles for L = 50 
at T* = 0.3 and N = 22 or fi* = -0.33904 (a), N = 23 or n* = -0.21127 (b), iV = 24 or 
= -0.02632 (c), N = 25 or = 0.22159 (d), N = 26 or = 1.08857(e). Repulsion to 
attraction ratio J* = 3 and attractive walls with hi = h*^ = —1 are considered for all the cases. 


IV. EQUATION OF STATE IN A SYSTEM CONFINED BY RIGID WALLS 


In this section we compute pressure for fixed number of confined particles Nq as a function 
of the distance between the conhning surfaces L. From these results we obtain the pressure 
as a function of density, p(p), for given Nq. For comparison we present p(p) calculated 
exactly in the GCE by the transfer matrix method described in Ref. 1^. In the GCE we 
consider hxed L and /i, and calculate p(/i) and p(/i) to obtain p(p) for given L. The shape 
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FIG. 2: The standard deviation of the number of particles ctn divided by the average number of 
particles {N) for L = 50, T* = 0.3, J* = 3 and attractive walls with = h’^ = —1. The red 
squares indicate the values of the chemical potential taken for the density profiles in Fig. [TJ 


(a) iji®3cocxxxxxxx:»«@a 

(b) BoooocaDoocccoooa 
ijiooocxxx»occc««®a 
ijiooocooooo(xx:»**a 

FIG. 3: Typical microscopic states present in the GS of the model with J* = 3, attractive walls 
with h\ = h*j^ = —1 and L = 15 in the CE with N = 5 (a), N = 6 (b) and N = 7 (c). For N = 5,7 
the GS is degenerate. 

of p{p) in the GCE depends on the commensurability between L and the period of the 
energetically favorable structure. We shall compare the results obtained in the CE for fixed 
Nq with the p{p) lines obtained in the GCE for six system sizes L. 

In Fig. |4]the results for the reduced pressure as a function of the average density for CE 
and GCE are presented. Note the discrepancy between the CE and GCE for p k, 0.55, where 
the periodically ordered clusters consisting of 3 particles are separated by 3 empty sites. In 
the GCE p[p) increases monotonically, although for p ^ 0.55 the slope is very large, and 
an inflection point is present. We should stress that in the CE Nq is fixed and the density 
changes because of the change of L. In contrast, in the GCE L is fixed, and the density 
changes because p, and as a result (iV), changes. In the GCE there are different branches of 
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p{p) for different L. One may interpret the nonmonotonic p{p) in the CE as a consequence 
of the jumps between the different branches of p{p) in the GCE for L and L — 1. 

In order to separate the effect of the fluctuation of the number of particles and the 
effect of the method by which the density changes, we compare the p{L) curves in the CE 
with Nq particles and in the GCE with (N) = Nq. In Fig. [5] the pressure is shown as a 
function of L for the CE with Nq = 21 and for GCE with {N) = 21. The GCE curve 


was obtained by Ending for each system size L the value of the chemical potentia 
that p{po) ~ Nq/L. We used the exact expression for density obtained in Ref. 


uq such 


19|. For 


such chemical potential the pressure was computed from the approximate formula /3p = 
— lnS(/io, T + 1, T) +lnS(/io, T, T), which is the Id lattice version of the standard expresion 
p = —{dfl/dV)^^T- We also present the density profiles for L = 35,36,37,38, where p 
changes rapidly in a nonmonotonic way. In the case of attractive surfaces the periodic 
structure where three occupied sites are separated by three empty sites is possible for L = 39, 
and corresponds to seven clusters. For L < 39 either the clusters are bigger, or the distances 
between them are smaller. In both cases the repulsion between the particles is present, 
and pressure increases. In Figl5]we can see six clusters for L < 37, and seven clusters for 
37 < L < 42. Note that the nonmonotonic dependence of p on L corresponds to the jump 
of the number of clusters. The unusual increase of pressure in the expanding system results 
from the transition to a larger number of smaller clusters. The clusters repel each other for 
L <39. Upon increase of the system size from L = 38 to L = 39 the separation between the 
clusters becomes large enough to put the clusters at the separations larger than the range 
of repulsion, and the pressure drops. 

It is interesting that although both the average densities and density profiles for L = 42 
in the two ensembles are the same, the pressure is different. The reason is that the pressure 
depends not on the values of the thermodynamic potentials at a given state, but on their 
change and as can be seen on panel (f) of Fig. [5l for L = 43 the profiles differ significantly. 

We conclude that the mechanical properties of a confined self-assembling system depend 
significantly on whether the system expands for given number of particles, or the separation 
between the system boundaries is fixed, and the number of particles decreases due to a 
change of the chemical potential. In both cases we can have the same change of density, 
but different changes of pressure. The unusual increase of pressure upon system expansion 
is found only in the case of fixed (average) number of particles, and is connected with a 
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significant structural reorganization. 




FIG. 4; Reduced pressure as a function of density for J* = 3, T* = 0.5 and attractive walls 
{h\ = h*^ = —1). Upper panel: Canonical Monte Carlo (CMC) simulation for the number of 
particles Nq = 25 and different system sizes L, obtained via virtual insertion method. Lower panel: 
GCE exact results for pressure vs. density for different system sizes, L = 50,51,52,53,54,55. 
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FIG. 5: Upper panel: reduced pressure as a function of the system size L for J* = 3 and T* = 0.5 
in the case of attractive walls {h^ = h*j^ = —1). Red curve with squares: Canonical Monte Carlo 
(CMC) simulation for 21 particles. Black curve with bullets: exact results in the GCE with the 
average number of particles equal to 21. Lower panels: (a-d) CE density profiles in systems with 
L = 35 (a), L = 36 (b), L = 37 (c) and L = 38. (indicated on the p{L) plot as blue triangles). 
Panels (e-f): comparison of the density profiles in CE (solid line) and GCE (dashed line) for L = 42 
(e) and L = 43 (f). 
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V. THE CASE OF ELASTIC BOUNDARIES 


In this section we assume that the separation between the system boundaries is not fixed, 
but can fiuctuate around L = Lq. The change of the wall separation is associated with the 
energy cost AU* = /c*(L —Lq)^. Here k* denotes the elastic constant in units l/a^. In Figl6] 
the average system size {L{N)) as a function of the number of particles N is presented. The 
confining surfaces are kept at the separation L by the spring that is at rest for Lq = 21. 
We assume attractive walls and T* = 0.5. In a system with rigid boundary conditions, 
attractive walls and L = 21, the periodic structure made of 12 particles is energetically 
favorable. Thus, for iV < 12 only for small values of the spring constant (e.g. k* = 0.1) 
significant deviations of the average system size from the reference value are present. For 
A > 12 the internal stress of the fluid competes with the elastic forces and (L) > Lq even 
for k* = 1. 

We can distinguish two limiting cases: (i) stiff spring, where the system size saturates 
and the particles become densely packed when N increases, and (ii) soft spring, where the 
average system size increases with increasing N, and the clusters are separated by empty 
sites. If k* is small enough, then the slope of {L{N)) increases when a new cluster made of 
3 particles is introduced to the system (see panels (a)-(f) of Fig. [6]). On the other hand, 
for larger values of k the elastic forces are stronger, and the system tends to modify the 
structure of the fluid rather then the system size. Between the two limiting cases there is 
an interesting region where the elastic and the solvation force are comparable and compete. 

In FigJT] we show {L{N)) for Lq = 21 and k* = 0.2 for several temperatures. At low 
temperatures three regimes with significantly different slopes of the lines {L{N)) can be 
distinguished. For A < 12 the slope of the line (L(A)) is small, because in this case the 
separation between the clusters ensures no repulsion between them. For 12 < A < 15 the 
system expands significantly upon addition of particles, because for A = 14 an additional 
cluster appears. In this region the average size of the clusters and the distance between them 
is 3. For A > 15 the slope is small again. Here elastic stress dominates and the clusters 
average size increases until the system becomes densely packed. Note that in this region and 
at low T, before the system gets filled with particles, we obtain an oscillatory dependence 
of (L(A)) on A, with the minima occuring when the number of particles is a multiple of 
3. Note also that for a given A, (L) increases with temperature except from 13 < A < 21, 
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FIG. 6: Average system size as a function of the number of particles for different spring constants 
k*. The walls are attractive, h\ = h\ = —1, J* = 3, T* = 0.5 and Lq = 21. Thick solid black 
lines, L = Lq and L = 2N — 3, correspond to rigid walls and to walls fully adapting to the optimal 
structure respectively. In the panels (a-f) the density profiles for the number of particls and the 
system size marked by the red squares along the black solid curve { k* = 0.1) are shown, (a) 
L = 23 and N = 12, (b) = 13 and L = 25, (c) L = 28 and N = 16, (d) L = 29 and N = 17, (e) 

L = 32 and N = 19, (f) L = 34 and N = 23. 


where for L = Lq the density is between the density of the periodic and the closely packed 
structures. We verihed that the anomalous contraction of the heated system is no longer 
observed at high T. Finally, note that there is some similarity of the shapes of the {L{N)) 
and p{p) lines (Figs 17] and 0]). 

The average wall separation and the average density prohle give insufficient information 
about the system behavior. In FigsU] and [U] we present histograms for the wall separation. 
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FIG. 7; Average system size as a function of the number of particles for different temperatures. 
Attracting walls assumed {h\ = h*^ = —1), J* = 3, spring constant k* = 0.2 and Lq = 21. The 
inset presents the average density {p) = N/{L) as a function the number of particles. 


Two cases can be distinguished - a single maximum in the probability of the appearance of 
the wall separation L, and a bistability with two maxima in this probability, separated by 
AL = 3. 



FIG. 8: The histograms for the wall separation in the case of elastic boundaries with the spring 
constant k* = 0.2, N = 16, Lq = 21, J* = 3 and attractive walls {h\ = /ij = —1) for T* = 0.2 (left 
panel), and T* = 0.5 (right panel). 

In order to understand the energetics associated with the bistability, let us consider the 
GS for Lq = 10 and N = 7. The microstates shown in Fig. flU] correspond to the same energy 
of the conhned system ([6]), when J* = 3 and k* = 2/9. Two different system sizes in the GS 
can occur when N is not a multiple of 3, and for L = Lq an intra-cluster repulsion is present. 
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FIG. 9: The histogram for the wall separation in the case of elastic boundaries with the spring 
constant k* = 0.06607. J* = 3, T* = 0.5, N = 17, Lq = 19, attractive walls {h\ = /i^ = — 1). 

The expansion is associated with a simultaneous increase of elastic energy of the walls, and 
decrease of the internal energy of the particles, when the separations between them are such 
that the repulsion is absent. Each microstate in Fig. [TUloccurs with the same probability, but 
because of the difference in the degeneracy for L = Lq and L = Lq + 3, the probability ratio 
for the two lengths is p{Lq + 3)/p{Lo) = 3. To estimate the spring constant leading in the 
above example to equal probability of Lq and Lq + 3 for low T*, we take into account only the 
microstates shown in Fig. [TUI and require that exp(—/9*(—5 + J*)) = 3exp(—/9*(—4 + 9k*)) 
(see (jU]) for L = Lq,Lq + 3). For T* = 0.2 we obtain k* ^ 2/9 + 0.0244136 in very good 
agreement with the results of simulations shown in Fig. [TUJ 
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FIG. 10: Upper panel: Two microstates in the GS of a system consisting of 7 particles with the 
repulsion to attraction ratio J* = 3, confined by attractive walls {h^ = = — 1) on a spring with 

the spring constant k* = 2/9 that is at rest for Lq = 10. In panel (a) L = Lq while in panel (b) 
L = Lq + 3. Note that the microstates symmetric to those also have the same energy, hence for 
L = 10 there are 2 different microstates with the same energy, and for L = 13 there are 6. These 
microstates correspond to the maxima of the probability shown in the lower panel, in which we 
present the histogram for the wall separation in the above system in the case of T* = 0.2 and 
elastic boundaries with the spring constant k* = 2/9 + 0.0244135. 
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VI. SUMMARY AND CONCLUDING REMARKS 


We have studied the effect of various constraints on colloidal self-assembly in thermody¬ 
namic states that correspond to self-assembly into small clusters separated by voids. In our 
model system small clusters with no intra-cluster repulsion yield a negative contribution to 
the internal energy, and do not interact with one another if the separation between them is 
larger than the range of repulsion. The positions, the size and the number of clusters can 
fluctuate and different deformations of the bulk structure are possible in confinement. In 
order to determine the role of constraints imposed on the number of particles and/or the size 
of the system, we have compared density profiles, equations of state, and effects of elastic 
boundary conditions. 

Structure: If the number of particles in the CE is too small for formation of the bulk 
structure, then one more cluster can be present in the GCE despite the same average number 
of particles (Fig[T]). Interestingly, the different number of clusters in the two ensembles is 
observed for thermodynamic states that do not correspond to the maximum of the fluctuation 
of the number of particles. This is because the fluctuations of the number of particles are 
coupled with the fluctuations of the size of the clusters. Even a small increase of the number 
of particles in the GCE together with the splitting of the clusters can lead to formation of 
a larger number of smaller clusters. The qualitative difference between the two ensembles 
occurs for quite large number of particles. This behavior is different than in simple fluids, 
where the difference between the two ensembles was observed for very small number of 


particles confined in very small pores 
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Equation of state: In the bulk the isotherms p{p) do not depend on the way in which 
the variation of density is attained. In the confined inhomogeneous systems it is no longer 
the case:/different curves p{p) are obtained when the size of the system changes with fixed 
number of particles, or when the number of particles changes at fixed system size. 

The shape of the p{p) curve obtained in the GCE with fixed L depends significantly on 
L, or more precisely on the commensurability of L and the period of the ordered structure 
(Fig H]). In the CE with fixed N we have obtained anomalous decrease of pressure for 
increasing density for small density intervals below and above the density of the equilibrium 
bulk structure (FigS]). Inspection of density profiles shows that the anomalous increase of 
pressure for increasing system size with fixed number of particles is accompanied by increased 
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number of clusters (the larger cluster splits) (FigE]). Recall that at short separations the 
clusters repel each other, and this leads to the increase of pressure. The pressure rapidly 
drops when L further increases and the clusters do not repel one another any more. In order 
to check if the anomalous dependence of pressure on density follows from the fixed number 
of particles or from the process by which the density varies, we computed the p{L) curve 
in the GCE with fixed average number of particles. We have obtained similar curves in the 
two ensembles except for large slits. In both ensembles the anomalous increase of pressure 
for increasing L is associated with the increasing number of clusters. For some large L the 
number of clusters increases with the system size in the GCE but not in the CE (Figj5]e,f). 
In this case a maximum in the p{L) curve is present only in the GCE. 

Bistability in elastic confinement: If the width of the slit can vary, then the system 
tends to equilibrate the competing solvation and elastic forces. We have found that the 
equilibrium size of the system is not always unique. In Ref. 19|] we observed a bistability in 
a system confined by elastic walls with permeable walls (fixed p). Two different system sizes 
can be equally probable: one with expanded and the other one with compressed boundaries. 
The size fluctuations are accompanied by an absorption or evaporation of a whole cluster. 
In the case of impervious walls (fixed N) a bistability exists too. In both cases the origin 
of the bistability is the change of the number of clusters, but the mechanisms which lead to 
the change are different. The number of clusters can fluctuate for fixed number of particles, 
when the state with a smaller number of bigger clusters and the state with a larger number 
of smaller clusters are equally probable. When the intra-cluster repulsion in large clusters 
competes with the elastic energy of stretched boundaries, the clusters can split and separate. 
The difference between the two equilibrium widths of the system is equal to the period of 
the bulk structure in the case of permeable walls, and to half the period of the bulk structure 
when the number of particles is constrained. Since in our Id system the energy barrier is 
~ ksT, spontaneous changes of the system size may be induced by thermal fluctuations. 

Our results show that different anomalies in the confined inhomogeneous systems occur 
when the release of some constraint or a change of the thermodynamic state leads to a change 
of the number of clusters. In particular thermodynamic states the structure and anomalies 
depend qualitatively on the ability of the system to interchange particles with its surroundings. 

Let us stress the difference between the confined simple fluids and the SALR systems. 
Packing effects of molecules or clusters play important role in both cases, especially for 
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the solvation force that exhibits oscillatory decay on the corresponding length scale in each 
case. However, the clusters can split or merge, unlike the molecules. This leads to qualitative 
differences between the ensembles, anomalies in the p{p) and the {L{N)) curves, and the 
bistability of the system size that in simple fluids are absent. 


The patterns emerging in the colloidal and amphiphilic self-assembly are very similar 
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. The clusters or layers composed of particles are distributed in space in a similar way 


as micelles or bilayers composed of amphiphilic molecules. Based on the similarity between 
the two types of self-assembly, we can expect that our results may also concern amphiphilic 
systems such as surfactants or lipids in water and block copolymers, at least on a qualitative 


level. Similarly, magnetic systems with competing interactions 
properties in confinement. 
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32| may have very similar 
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Appendix A: The expressions for H* and Y1's(m) 


The H* contains the pair interactions between the particles at the sites of the M-th box, 
as well as the interactions between the particles at the sites labeled 3M-I-1 and 3M-I-2 (if 
they exist for given L). 


H*[S{M)] 


-(EloP(3M-^)p(3M-^-l))-p*(EloP(3M-^)), ifj = 0 

-(ELoP(3^+1 “ i)p{3M-i)) + J*p{3M-2)p{3M+1) if j = l 

^ -/^*(EloP(3M + l-*)), 

-(E?=oP(3^+2 - i)K3M+ 1 - i))+ if j = 2 

, ^*(EIoP(3M-2 + t)pi3M+l + t)) - p*(EtoP(3M + 1 - z)). 
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Whereas Y1's{m) denotes 


S{M) 


= 


J2siM) if j = 0 

J2s{M) X]p(3M+l) if J = 1 

J2s{M) X]p(3M+l) X]p(3M+2) if i = 2 


Appendix B: Microstates obtained by the virtual expansion of the system 


We consider N indistinguishable particles and L lattice sites. Each site can be empty or 
occupied by one particle, thus there are (^) distinguishable microstates. We will show that 
the virtual expansion procedure of building configurations of the system with L + 1 sites by 
inserting an empty site at a random position in a system with L sites is not biased by the 
insertion procedure. Let us consider two sets of particle configurations. The elements of the 
first set are the microstates of a system of size L with N occupied sites with a distinguished 
position in which a new site can be inserted. Since the new site can be inserted in L + 1 
places, the first set contains (L + 1)(^) different elements. In the second set each of the 
elements corresponds to one of the microstates of a system of L + 1 sites with N occupied 
sites, and with one of its L + 1 — iV empty sites marked as removable. The number of 
elements of the second set is (L + 1 — iV) ■ Because {L + 1) = (L + l — N) , the 

two sets are equinumerous. We define a one-to-one correspondence between the elements of 
the two sets by identifying the location of the insertion point in an element of the first set 
with the location of the removable site in the element of the second set, and by requiring 
that the same sites are occupied. Note that each microstate of the system with L -|- 1 sites 
can be obtained in L -|- 1 — A ways from the elements of the second set by removing the 
mark “removable”. Because of the one-to-one correspondence between the elements of the 
two sets, each microstate of the system with L -|- 1 sites is obtained L + 1 — N times by the 
above procedure. 

It follows that the proposed procedures of building configurations of the system with 
L -|-1 sites by inserting an empty site at a random position on the configurations of a system 
with L sites are not biased by the insertion procedure. The same lack of bias applies in the 
reverse procedure. 

For an illustration let us consider L = 2 and iV = 1. There are 2 microstates, [1,0] and 
[0,1]. After insertion of an empty site in 3 possible places, we obtain from the first microstate 
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[0,1,0], [1,0,0], [1,0,0], and from the second microstate [0,0,1],[0,0,1], [0,1,0]. One can easily 
see that after this procedure we obtained each microstate in the system of size L + 1 = 3 
containing N = 1 particle L + 1 — iV = 2 times. 
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